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Abstract 

The (Michelson) Sivashinsky equation of premixed flames is studied in a rectangular domain in 
two dimensions. A huge number of 2D stationary solutions are trivially obtained by addition of two 
ID solutions. With Neumann boundary conditions, it is shown numerically that adding two stable 
ID solutions leads to a 2D stable solution. This type of solution is shown to play an important role 
in the dynamics of the equation with additive noise. 
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I. INTRODUCTION 

The Sivashinsky equation jl| (or Michelson Sivashinsky equation depending on the au- 
thors) is a non linear equation which describes the time evolution of premixed flames. Be- 
cause of the jump of temperature (and thus of density) across the flame, a plane flame front 
is submitted to a hydrodynamic instability called the Darrieus-Landau instability. The con- 
servation of normal mass flux and tangential velocity across the front leads to a deflection 
of streamlines which is the main cause of this instability. A more detailed description of 
this instability can be found in the book of Williams [2|| (see also, in the approximation 
of potential flow in the burnt gases, the elementary electrostatic explanation in Q|, where 
essentially the flame is described as a surface with a uniform charge density). At small 
scales, the instability is damped by diffusive effects: the local front propagation velocity is 
modified by a term proportional to curvature, the coefficient ahead of the curvature term 
is called the Markstein length. A geometrical non linear term, which limits ultimately the 
growth of the instability, is caused by the normal propagation of the flame. The Sivashinsky 
equation, obtained as a development in powers of a gas expansion parameter, i.e. for a small 
jump of temperature, or equivalently for a flow almost potential in the burnt gases, repre- 
sents a balance between the evolution due to these three effects, Darrieus Landau instability, 
stabilization by curvature and normal propagation of the flame. 

The qualitative agreement between the Sivashinsky equation and direct numerical simu- 
lations, generally performed with periodic boundary conditions has been excellent even with 



large gas expansion Q], and also when gravity is included 

paper of the field [f| (following [zj], where the pole decomposition was introduced) that the 



5|. It has been shown in a classic 



ID solution of the Sivashinsky equation in the absence of noise was attracted for large times 
toward stationary solutions, with poles aligned in the complex plane, called coalescent solu- 

no 

tions. It was shown analytically in [8||9fl that each solution, with a given number of poles, 
is linearly stable in a given interval for the control parameter (either the domain width or 
more often the Markstein length with a domain width fixed to 27r). 

In a recent paper |l(3] (hereafter called I) the present author has been interested in the 
behavior of the Sivashinsky equation in ID, but with Neumann boundary conditions (zero 
slope at each end of the domain), a situation which, although more realistic than periodic 
boundary conditions, had not attracted much interest over the years. Actually, periodic 
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boundary conditions lead to a symmetry which is not present in the case of a flame in a 
tube, i.e. every lateral translation of a given solution is also a solution. Presented in a 
different way, a perturbation on the flame can grow, reach the cusp (the very curved part 
of the front, pointing toward the burnt gases) and then decay, but after having caused a 
global translation of the original solution. This is not possible with Neumann boundary 
conditions, but it was supposed that this difference with periodic boundary conditions was 
unimportant. The surprise was however that stable stationary solutions in the Neumann 
case involved a number of solutions with two cusps (and the corresponding poles) at each 
end of the domain, called bicoalescent solutions. This type of solution of the Sivashinsky 
equation was already introduced in [11], although this last article did not obtain those which 
are stable with Neumann boundary conditions^ The author would like to mention here two 
articles which he did not cite in I , namely 



12 ], where some bicoalescent solutions with 
Neumann boundary conditions were first obtained, and {3] where bicoalescent solutions 
were obtained in direct numerical simulations. In this last paper, one solution was not 
completely stationary, because of the effect of noise, but another solution was actually 
almost stationary. Of course the computer time needed for such a simulation is probably 
one hundred times more than the equivalent Sivashinsky equation simulation, with all sorts 
of possible sources of noise, so obtaining really stationary bicoalescent solutions in this case 
is a challenging task. 

Coming back to I, we can summarize the ID results of this paper in the following way: 

1. Bicoalescent solutions were obtained, stable with Neumann boundary conditions. Sim- 
ulations performed without noise tend to these solutions. 

2. The new solutions led to a bifurcation diagram with a large number of stationary 
solutions, where particularly the number of solutions multiply when the Markstein 
length, presented above, which controls the stabilizing influence of the curvature term, 
decreases 

3. The bicoalescent solutions play a major role in the dynamics of the equation with 
additive noise. In the case of moderate white noise, the dynamics is controlled by 
jumps between different bicoalescent solutions. 

In the present paper, we shall be interested in the Sivashinsky equation with Neumann 
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boundary conditions, but in two dimensions in a rectangular domain. Another nice property 
of the equation (apart from the pole decomposition in ID) is that 2D solutions can be 
formed by the simple addition of two ID solutions, one for each coordinate jlj]]. The exact 
counterpart of I will be obtained: 

1. Sums of two bicoalescent solutions are stable in 2D with Neumann boundary condi- 
tions. The time evolution of the equation without noise tends toward these solutions 

2. With sums of a large number of ID solutions, a really huge number of 2D solutions 
can be obtained. 

3. The sums of bicoalescent solutions play also a major role in the dynamics in two 
dimensions in the presence of noise. 



II. SOLUTIONS IN ELONGATED DOMAINS 

The Sivashinsky equation in one dimension can be written as 

<Pt + \<t>l = H** + I {4>) (1) 

where (ft (x) is the vertical position of the front. The Landau operator / (0) corresponds 
to a multiplication by \k\ in Fourier space, where k is the wavevector, and physically to the 
destabilizing influence of gas expansion on the flame front (known as the Darrieus-Landau 
instability, and described in the introduction), v is the only parameter of the equation (the 
Markstein length) and controls the stabilizing influence of curvature. The linear dispersion 
relation giving the growth rate a versus the wavevector is, including the two effects: 

o = \k\ - uk 2 (2) 

As usual with Sivashinsky-type equations, the only non linear term added to the equa- 
tion is \(j> 2 x . In the flame front case, this term is purely geometrical : the flame propagates 
in the direction of its normal, a projection on the vertical (y) direction gives the factor 
cos (8) = l/\Jl + 4>l, where 9 is the angle between the normal and the vertical direction, 
then a development valid for small slopes of the front leads to the term \4>\- The Sivashinsky 
equation is typically solved numerically on [0, 2tc] with periodic boundary conditions. In I 
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it has also been solved on [0, 2tt] with only symmetric modes, which corresponds to homo- 
geneous Neumann boundary conditions on [0, ir] (zero slope on both ends of the domain). 
The two dimensional version of the Sivashinsky equation is 



<Pt + \w) 2 = uA<p + (3) 



where the Landau operator / ((f)) corresponds now to a multiplication by ^Jk x + ky in 
Fourier space, k x and k y being the wavevectors in the x and y directions. All dynamical 
calculations, are performed by Fourier pseudo-spectral methods (i.e. the non linear term is 
calculated in physical space and not by a convolution product in Fourier space) . The method 
used is first order in time and semi-implicit (implicit on the linear terms of the equation, 
explicit on \4> 2 X ). No particular treatment of aliasing errors is used. The 2D Sivashinsky 
equation is solved in [0, 2n] * [0, 26] with only symmetric modes, which corresponds to ho- 
mogeneous Neumann boundary conditions in the rectangular domain [0, ir] * [0, b] . 

Pole solutions (jf|) of the ID Sivashinsky equation are solutions of the form: 

. 2 ,|j ln ( sin (i^)) +b ( sb (^W))} (4) 

where N is the number of poles z n (t) in the complex plane. Actually the poles appear in 
complex conjugate pairs, and the asterisk in Equation [4] denotes the complex conjugate. In 
all the paper, the number of poles will also mean number of poles with a positive imaginary 
part. The pole decomposition transforms the solution of the Sivashinsky equation into the 
solution of a dynamical system for the locations of the poles. In the case of stationary 
solutions, the locations of the poles are obtained by solving a non linear system: 

-v ^2 cotf^— — - J - isgn [Im (z n )\ = n=l,---,N (5) 

where Im (z n ) denotes the imaginary part and sgn is the signum function. This non linear 
system is solved by a variant of Newton method. 

With periodic boundary condition, the usual result is that in the window 2n — 1< \jv< 
2n + l,n = l,2,-- - there exists n different monocoalescent stationary solutions (all the poles 
have the same real part), with 1 to n poles, and the solution with the maximum number of 
poles n is asymptotically stable. For a particular value of 1/z/, the number n(u) such that 
2n — 1< \jv< 2n + 1 is called the optimal number of poles. 



With Neumann boundary conditions, in each of the intervals [2n — 1, 2n + 1] of the pa- 
rameter 1/u, not only one asymptotically stable solution, but n + 1, of the form (l,n — I) 
with I — 0, 1,- • • ,n where I poles coalesce at x = and l — n coalesce at x = %, were obtained 
in I. (The bicoalescent type of solutions have been recently introduced in fUj). In Figured] 
is shown a bifurcation diagram with all the possible stable stationary solutions (plotted only, 
contrary to I, in the domain where they are stable) versus 1/u. What is actually plotted is 
the amplitude A0 (maximum minus minimum of 0) versus 1/u. As can be seen, when the 
optimal number of poles increases with 1/u, the number of stable stationary bicoalescent 
solutions is also increasing. The stability of these solutions is not proved analytically, nor 
by a numerical study of the linearized problem, we use only numerical simulations of the 
Sivashinsky equation, with the different bicoalescent solutions plus some small perturbations 
as initial conditions, and the solution returns toward the unperturbed solution. 

In a square domain [0, 2ir] * [0, 2n], it has been remarked in jlj]] that if (f>i(x) and faix) 
are solutions of the ID Sivashinsky equation (pQ), then <f>i(x) + 02 (y) (we use here 0i © 02 
as a notation for this sum, whose amplitude is the sum of the amplitudes of 0i and 02) is 
a solution of j3]) in two dimensions, and that the stationary solution obtained numerically 



in this case for periodic boundary conditions [15| is simply a sum of two monocoalescent 
ID solutions. Let us note that, if it is absolutely obvious that sums are solutions of the 2D 
equation, the stability of these solutions has never been proved analytically, and can only 
be inferred from a small number of numerical simulations. 

In the case of rectangular domains [0, 2ti] * [0, 26], sums are also solutions of the equation, 
with 02 now solution of the ID Sivashinsky equation with parameter 1/u in domain [0,26], 
which can be obtained by an appropriate rescaling from the solution in [0, 2ir] with parameter 
1/u, = (l/u)(b/n). 

A particularly simple case is the limit where 6 is very small, where the only solution with 
parameter 1/u in domain [0,26] is simply the flat (0) solution 2 = 0. As a sum of the 
previously described bicoalescent solutions in [0, 2n] added to the flat solution in the other 
direction, we have simply a way to observe the bicoalescent solutions in two dimensions. We 
have observed numerically (not shown here, the behavior is very similar to the ID case) for 
Neumann boundary conditions, that these sums (l,n — I) © (0) are stable . As an example, 
for 1/u = 10 and 6 = n/10 we show in Figure [2] a perspective view of the three different 
stationary bicoalescent solutions (5, 0) © (0), (4, 1) © (0), (3, 2) © (0) (from top to bottom). 
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In all the figures, the whole domain [0, 2ti} * [0, 26] is plotted, the solution with Neumann 
boundary conditions corresponds only to one fourth of the domain [0, it] * [0, b] . We have 
found it clearer to show the whole domain (contrary to I), because some solutions are very 
difficult to distinguish if plotted in [0, 7r] * [0, b}. Although these solutions are very sensitive 
to noise (although less than the pure ID solutions) it could be possible to observe in direct 
numerical simulations and experimentally the solutions with the lower amplitude, which are 
the least sensitive to noise. In experiments, the solutions should also survive heat losses 
(important in narrow channels) and not be too much perturbed by gravity (i.e. have a large 
enough Froude number) in order to be observed . 

III. SOLUTIONS IN SQUARE DOMAINS 

We now turn to stationary solutions in square domains [0, 27r] * [0, 2ti] with Neumann 
boundary conditions. Sums of bicoalescent solutions produce also in this case stable sta- 
tionary solutions. The purpose of this section is to give details on the consequences of this 
simple addition property. We show first the different types of solutions obtained by addition 
of stable bicoalescent solutions in ID. We insist on the fact that these solutions are linearly 
stable and give a specific example of the time evolution of one such solution with some small 
perturbations. Finally two bifurcation diagrams are provided, one is the 2D equivalent of 
Figure Q] with the stable solutions plotted only in their stable domain. The second contains 
all the solutions obtained by addition of all the branches found in ID in I, and as the reader 
will really huge number of branches are created in this way. 

In Figures [3] and [Hare shown the six stable solutions obtained from the three ID solutions 
of Figure [2] for 1/v — 10. In Figure [3] can be seen (in perspective view, for the whole domain 
[0,27r]* [0,27r]), from top to bottom the (5,0)0(5,0), (4, l)©(5,0)and (3, 2)©(5,0) solutions. 
In Figure [4] can be seen the three remaining solutions (4,1) © (4,1), (3,2) © (4,1) and 
(3,2) © (3,2). All these solutions are found to be linearly stable, although all the solutions 
of Figure[3] ((5, 0)© something) are extremely sensitive to noise. It must be pointed out that 
most of these solutions would have been almost impossible to find from a time integration 
of the 2D Sivashinsky equation (Equation ([3])) because of this sensitivity to noise, and it is 
likely that obtaining them from a steady version of j3j) would have been very difficult too. 

In Figure [5], we have an example showing the stability of the (3, 2) © (4, 1) solution. We 
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start from this solution and add an additive white noise to Equation j3]) when the time is 
below 0.5. This white noise is gaussian, of deviation one, and we multiply it by an amplitude 
a = 0.001. It can be seen that after the noise is stopped, the solution tends exponentially 
back toward the (3, 2) © (4, 1) solution. Similar figures would be obtained with the other 
solutions of Figures [3] and SI except that higher amplitude solutions would need an even 
lower noise in order not to jump immediately toward a lower amplitude solution. 

In Figure [6] is shown the strict 2D equivalent of Figure [TJ the bifurcation diagram showing 
the amplitude versus 1/v for all the solutions which are linearly stable, only plotted in their 
domain of stability. For 1/v < 3 there is only one possibility (1, 0) © (1, 0). For 3 < 1/V < 5 
we have three branches (from higher to lower amplitudes) (2,0) © (2,0) (1, 1) © (2,0) and 

(1.1) © (1,1). For the value 1/V = 10 we have the six solutions of Figures [3] and 31 that is 
from higher to lower amplitudes the (5, 0) © (5, 0), (4, 1) © (5, 0) (3, 2) © (5, 0) (4, 1) © (4, 1), 

(3.2) © (4,1) and (3,2) © (3,2) solutions. Higher values of 1/V would correspond to an 
increasing number of stable stationary solutions. 



Naturally, neither Figured] (in ID) or 
solutions. In ID Guidi and Marchetti 



jffure M (in 2D) contain all the possible stationary 



ll| have introduced the concept of interpolating 



solutions, which are unstable solutions connecting different branches of stable solutions in 
the previous bifurcation diagrams. In I, the present author has shown that this leads to a 
complex network of solutions, which was called web of stationary solutions. But now in two 
dimensions, we have the possibility, when two branches 0i and 02 exist for a parameter 1/V 
to create the 2D branch <pi © 2 . This construction leads to a bifurcation diagram (with as 
before 1/v < 14, i.e. not very large flames) with a truly huge number of different stationary 
solutions (several thousands of branches). The comparison with Figure [6] shows that most 
of these solutions are linearly unstable. 

The author would like to insist here on different points. First, it is only possible to obtain 
such an incredible number of stationary solutions because of two properties of the Sivashinsky 
equation, the pole decomposition, which transforms the search of stationary solutions in one 
dimension in a 0D problem, then the possibility to add ID solutions in order to get 2D 
rectangular solutions of the Sivashinsky equation. In the Kuramoto-Sivashinsky equation 
case (a non linear equation with a different growth rate but the same non linear term) 
the pole decomposition is not available, but nevertheless a lot of ID stationary solutions 
have been obtained |lf3|. The Kuramoto-Sivashinsky equation shares with the Sivashinsky 
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equation the possibility to create 2D solutions by adding two ID solutions, so actually in this 
case we have also a very large number of branches. These rectangular solutions are not as 
physically relevant in the Kuramoto-Sivashinsky equation case. Contrary to the Sivashinsky 
equation, where stable stationary solutions are basically as large as possible and are thus 
rectangular in a rectangular domain, it seems likely that in the Kuramoto-Sivashinsky case, 
the most interesting solutions would have an hexagonal symmetry (hexagonal cells are also 
observed for the Sivashinsky equation with stabilizing gravity 13]). Stationary solutions 
of the Sivashinsky equation with hexagonal symmetry should exist too, and the author 
conjectures that the order of magnitude of the number of solutions with hexagonal symmetry 
should be approximately the same as those with rectangular symmetry. Apparently there 
is no trivial way to construct hexagonal solutions, so unfortunately, until some progress is 
made, obtaining the hexagonal equivalent of Figure [7] is almost impossible. We have here 
an example emphasizing the fact that as the smoothing effect (viscosity, curvature, surface 
tension ...) decreases, we are not able to generate correctly all the simple solutions of a 
given set of partial differential equations (Sivashinsky and Kuramoto-Sivashinsky equations, 
Navier Stokes ...) even with the aid of computers. 



IV. EVOLUTION WITH NOISE 



In the previous section we have shown numerically that the sums of linearly stable ID 
bicoalescent solutions lead to linearly stable 2D solutions. However, even a linearly stable 
solution could have a very small basin of attraction. So in this section, we study the effect 
of noise on the solutions of the Sivashinsky equation in a square domain, with Neumann 
boundary conditions. The important solutions will be the solutions that are reasonably 
resistant to the applied noise. 

This noise used here is simply an additive noise, added to the right-hand side of Equation 
([3j). We choose the simplest possible noise, a white noise (in space and time), which is 
gaussian, has deviation one and is multiplied by an amplitude a. But contrary to Figure 
El this noise will be applied at each time step. We use in all the simulations presented the 
same parameter \/v = 10, the stationary solutions corresponding to this parameter have 
been presented in the previous section. We recall that in I, for the one dimensional version 
of the Sivashinsky equation with moderate noise, the evolution was analysed in terms of 
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jumps between the available bicoalescent stationary solutions. We would like to show here 
that in 2D, the sums of bicoalescent solutions also play an important role in the dynamics. 

In FigureEl starting from an initial condition which is the (4, l)ffi(4, 1) stationary solution, 
is plotted the amplitude of the solution versus time, for a noise amplitude a = 0.01, with also 
straight lines corresponding to the amplitudes of the lowest amplitude stationary solutions, 
i.e. those of Figure [H The stationary solutions with higher amplitudes (those of Figure [3D 
apparently are too sensitive to noise to play any role in the dynamics. It is seen in Figure 
[8] that because of the noise, the solution departs quickly from the (4, 1) © (4, 1) solution, 
and that it seems that, during the time evolution, the solution is close (apart from some 
violent peaks in the amplitude) to the (3, 2) © (4, 1) solution for some time, then finally the 
amplitude decreases again to be near that of the (3, 2) © (3, 2) solution. 

In order to prove that the solution is indeed close to the previously mentioned solutions, 
because after all very different solutions could have similar amplitudes, we plot in Figure [9] 
for the same simulation, what we have called the distance between the solution at a given 
time and the sums of bicoalescent solutions, which is simply the LI norm of the difference 
between both solutions. The spatial mean value of all solutions is adjusted here to have 
the same value. Normally, it is necessary to measure the distance between the solution and 
all symmetries of a given sum of bicoalescent solutions (i.e. you can interchange the poles 
at and 7r in the x and y directions) but for the low amplitude a = 0.01, it has not been 
necessary, and we plot only the distance from the relevant solutions. 

As we start from the (4, 1) © (4, 1) solution, the distance to this solution is zero initially, 
and we can see that, although the amplitude seems to indicate that at some time, one is again 
close to this solution, this is not the case. On the contrary, the solution returns regularly 
close to the (3, 2) © (4, 1) solution for times lower than 110, then there is a transition toward 
something close to the (3, 2) © (3, 2) solution, the solution departs slightly from this last 
solution for some time, possibly toward a linearly unstable stationary solution, and returns 
toward it at the end of the simulation. As Figure [6] remotely looks like energy levels in 
atomic physics, one could be tempted to interpret the evolution of the two previous figures 
with a small noise (apparently in 2D the solution is less sensitive to a given amplitude of the 
white noise compared to ID simulations) as a sort of deexcitation from the high amplitude 
level (4,1) ©(4,1) toward first (3, 2) ©(4,1), then toward the fundamental level (3, 2) ©(3, 2). 
Indeed, between the sums of bicoalescent solutions, if all are linearly stable, the solutions 
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with the lower amplitude seems to be more resistant to the action of noise. 

To better understand the effect of noise, we present now a simulation with a larger 
noise amplitude a = 0.1 , ten times larger than the previous case (we recall that this noise 
amplitude should be compared to the laminar flame velocity, which is normalized to 1 in 
this paper). In Figure d0] is plotted the amplitude versus time, with as before straight lines 
with the amplitude of the important sums of bicoalescent solutions. The initial condition is 
also the (4, 1) © (4, 1) solution. Apparently this last solution is too sensitive to noise to play 
a meaningful role in the dynamics, although it happens that some peaks in the amplitude 
could involve solutions not too far from this initial solution. As the distance to this solution 
is never really small, even in the peaks, we shall not comment further on this solution. On 
the other hand, it seems that a lot of time is spent with an amplitude close to that of the 
(3, 2) © (3, 2) solution (which we have called previously the fundamental level), and perhaps 
some time with an amplitude close to the (3, 2) © (4, 1) solution (the first excited level). 

In order to see what is really occurring, we now turn to figures of the distance (defined 
above) to these two solutions versus time (for the same simulation of Figure fTOll. However, 
for a higher amplitude, we have to include the four different symmetries of these solutions 
in the analysis (i.e. for instance (3,2) © (4,1) (3,2) © (1,4) (2,3) © (4,1) (2,3) © (1,4) ) 
. In Figure [TT] is shown the distance to the four symmetries of the (3, 2) © (3, 2) solution. 
The distance to one of the four symmetries is indeed often small (but not very small for this 
value of the noise) during the time evolution. Then because of the noise, perturbations are 
created that lead the amplitude to increase as the perturbation is convected toward one of 
the cusps, and the solution often comes back toward another symmetry of the fundamental 
level. 

In Figure [12] is shown the distance to the four symmetries of the (3, 2) © (4, 1) solution 
(the first excited level) (always for the same simulation). It is seen that the solution is 
only reasonably close to this type of solution at times close to 50. At other times, minima 
of the distance are not very small and the solution is often closer to the (3, 2) © (3, 2) 
solution. In I, we have presented the evolution of the ID Sivashinsky equation with a 
moderate additive noise as a series of jumps between bicoalescent solutions. In 2D the 
situation is relatively similar, with the sums of bicoalescent solutions playing the same 
role. However, the noise amplitude necessary to cause jumps seems much higher in 2D, and 
practically speaking during the previous simulation, only the fundamental (the solution with 
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the lowest amplitude) and first excited levels were obtained. It should also be noted that the 
degenerescence (the four possible symmetries) of the fundamental level is probably important 
in the evolution (for instance for ljv = 12 the fundamental level would be (3,3) © (3,3), 
which does not lead to other solutions by symmetry, so that it should be less likely to obtain 
the fundamental level in this case). 

Before closing this section, let us insist on the fact that, if the solution regularly returns 
toward sums of bicoalescent solutions, the fronts we obtain are not sums for each time. Figure 
[T3l where is plotted a front of the previous simulation for time 120.1 (just before a peak of 
the amplitude in Figure \T0§ , should be a clear example of this property. In this Figure, 
the whole domain [0, 2n] * [0, 27r] is plotted as before, but this time as a grayscale figure, 
white corresponding to the minimum of (p, black to the maximum. Essentially, an oblique 
perturbation has grown on a front that was previously a sum. This oblique perturbation 
moves toward each corner of Figure HU and the amplitude peak corresponds to the moment 
where the perturbation reaches the corner. Then the solution is attracted again toward a 
sum of bicoalescent solutions. 

To summarize this section on the effect of noise, the fact that all the sums of bicoalescent 
solutions with the optimal number of poles are linearly stable does not prove that they can 
be practically observed. On the contrary, the solutions with the larger amplitude have a 
basin of attraction so small that they can almost never be seen. We have introduced an 
analogy with atomic physics by calling the bicoalescent solution with the lowest amplitude 
the fundamental level, other solutions the excited levels. In the examples shown, only the 
fundamental and first excited levels (and their symmetries) were obtained during the time 
evolution of the Sivashinsky equation excited by an additive noise. We recall that in I, it was 
shown in ID that the evolution with noise was completely different with periodic boundary 
conditions, where only the largest amplitude monocoalescent solution was linearly stable 
(even if extremely sensitive to noise). In this case, the solution regularly returns close 
to the highest amplitude solution. With Neumann boundary conditions, this is just the 
opposite, the solution prefers to be close to the lowest amplitude, almost symmetric, sum of 
bicoalescent solutions. 
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V. CONCLUSION 



In this paper, we have used the possibility to create two dimensional rectangular sta- 
tionary solutions from the addition of two ID stationary solutions in order to generate a 
huge number of stationary solutions of the Sivashinsky equation. With Neumann boundary 
conditions, the addition of two stable ID bicoalescent solutions leads to stable 2D solutions, 
which also play a role in the dynamics when an additive noise is added to the equation. How- 
ever, with noise, only the sums of bicoalescent solutions with the lowest amplitude (which 
are less sensitive to noise) have a reasonable chance to be observed. More precisely, jumps 
between different symmetries of the lowest amplitude sum, or between the two sums with 
the lower amplitude, are obtained in the simulations. Although we have used a white noise 
in this paper, experiments, submitted to a residual turbulence, should behave in a similar 
way. In order to have a large enough Froude number for gravity effects to be negligible, 
flames with a sufficiently large laminar flame velocity would have to be chosen. 
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Figure 1: Stable stationary solutions in ID: amplitude Acj) vs 1/u. All the different branches are 
only plotted for the values of 1/u where they are stable. A notation like (3,2) means that 3 poles 
are located at x = 0, and 2 poles at x = ir 
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Figure 2: Perspective view of the (from top to bottom) (5,0)ffi(0), (4, l)ffi(O), (3,2)©(0) stationary 
solutions for 1/v = 10 and b = n/10. The solution is plotted in the interval [0, 2ir] * [0, 26] because 
it is easier to visualize. Actually, Neumann boundary conditions are satisfied in [0, tt] * [0, b] (one 
fourth of the domain shown). 
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Figure 3: Perspective view of the (from top to bottom) (5,0) (5,0), (5,0) © (4, 1), (5,0) © (3,2) 
stationary solutions for 1/v = 10. The solution is plotted in the interval [0, 2ir] * [0, 2tt] because 
it is easier to visualize. Actually, Neumann boundary conditions are satisfied in [0, ir] * [0, n] (one 
fourth of the domain shown). 
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Figure 4: Perspective view of the (from top to bottom) (4, 1) (4, 1), (4, 1) (3, 2), (3, 2) © (3, 2) 
stationary solutions for 1/v = 10. The solution is plotted in the interval [0, 2ir] * [0, 2tt] because 
it is easier to visualize. Actually, Neumann boundary conditions are satisfied in [0, ir] * [0, n] (one 
fourth of the domain shown). 
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Figure 5: Amplitude vs time for \jv = 10 , starting from a (4,1) © (3,2) solution. A gaussian 
white noise (amplitude a = 0.001 ) is imposed on this solution when time is smaller than 0.5. The 
solution returns exponentially toward the initial solution. 
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Figure 6: Stable stationary solutions in 2D for a square domain : amplitude Acft vs 1/v. All the 
different branches are only plotted for the values of 1/V where they are stable. The 2D linearly 
stable solutions are obtained by addition of the corresponding ID linearly stable solutions of Figure 

m 
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Figure 7: Stationary solutions in 2D for a square domain : amplitude A(p vs 1/v (figure with all 
the solutions obtained by addition of ID stationary solutions). When two ID branches found in I 
coexist for a certain value of 1/v, a 2D branch is created, whose amplitude is the sum of the ID 
amplitudes. 
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Figure 8: Amplitude vs time for 1/u = 10 and a = 0.01 (low noise amplitude). Deexcitation from 
the (4, 1) © (4, 1) solution toward the (3, 2) © (3, 2) solution. This diagram suggests that the solution 
is first close to the (3, 2) © (4, 1) solution, then from the (3, 2) © (3, 2) solution, i.e. that the solution 
with the lowest amplitude is the most noise resistant. 
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Figure 9: Distance to the main stationary solutions vs time for 1/u = 10 and a = 0.01. A distance 
is a norm of the difference between the solution at a given time and the stationary solution. This 
diagram makes it possible to verify if a solution at a given time is indeed close to a stationary 
solution. 
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Figure 10: Amplitude vs time for 1/v = 10 and a = 0.1 (moderate noise amplitude). This figure 
suggests that the solution is often close to the (3, 2) © (3, 2) solution. It will be shown in the 
following figures that it is close to the (3, 2) © (4, 1) solution only for times around 50. 
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Figure 11: Distance to the different symmetries of the (3,2) © (3,2) solution vs time for 1/u = 10 
and a = 0.1. The noise is sufficiently large to induce transitions between the different symmetries 
of this fundamental level. 
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Figure 12: Distance to the different symmetries of the (3, 2) © (4, 1) solution (first excited level) vs 
time for 1/v = 10 and a = 0.1. The solution is only close to one symmetry of the (3,2) © (4,1) 
solution for times around 50 (after apparently a transition from a (3,2) © (3,2) solution). 
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Figure 13: Solution at time 120.1 for 1/v = 10 and a = 0.1, plotted as a grayscale figure (white: 
minimum of </>, black: maximum of <fi). The solution is plotted in the interval [0, 2ir] * [0, 2ir]. 
Presence of an oblique perturbation which has grown on a sum of bicoalescent solutions. This 
perturbation will reach one corner in the figure, be damped, and the solution will again be close to 
a sum of bicoalescent solutions. 
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